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ABSTRACT 

Aims. We present kinematic and photometric evidence for an accretion event in the halo of the cD galaxy M87 in the last Gyr. 
Methods. Using velocities for ~ 300 planetary nebulas (PNs) in the M87 halo, we identify a chevron-like substructure in the PN 
phase-space. We implement a probabilistic Gaussian mixture model to identify PNs that belong to the chevron. Prom analysis of deep 
V-band images of M87, we find that the region with the highest density of chevron PNs is a crown-shaped substructure in the light. 
Results. We assign a total of NpN sui, = 54 to the substructure, which extends over ~50 kpc along the major axis where we also observe 
radial variations of the ellipticity profile and a colour gradient. The substructure has highest surface brightness in a 20kpc x 60kpc 
region around 70 kpc in radius. In this region, it causes an increase in surface brightness by >60%. The accretion event is consistent 
with a progenitor galaxy with a V-band luminosity of L = 2.8 ± 1.0 x IO^Lq.v, a colour of (B - V) = 0.76 ± 0.05, and a stellar mass 
ofM = 6.4±2.3 X IO^Mq. 

Conclusions. The accretion of this progenitor galaxy has caused an important modification of the outer halo of M87 in the last Gyr. By 
itself it is strong evidence that the galaxy’s cD halo is growing through the accretion of smaller galaxies as predicted by hierarchical 
galaxy evolution models. 
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1. Introduction 


According to the current theory of hierarchical formation and 
evolution of structures, accretion events are believed to have 
an important role in the cosmological build up of stellar halos 
in elliptical galaxies dPe Lucia & Blaizolll200'^. responsible fo r 
their growth at relatively low redshifts (z < 2: IOser et al.llTOloh . 
In dense environments accretion is even more dramatic, such 
that close to the dynamical centre of the cluster, central cluster 
galaxies are expected to have the maiori ty of their stars accreted 
dLaporte et al.ll201^ ICooper et al.ll2014l) . 

On the observational front, the outer stellar envelopes are 
observed to increase in m ass by a factor of ~ 4 since z - 2 
(Ivan Dokkum et all l2010l) . The presence of an accreted com¬ 
ponent is usually identified as an excess o f light over the ex¬ 
trapolation of the g a Lxy’s inner profile (IZibetti et al.l l2005l 
iGonzalez et al.l 1200'^ D’Souza et al. or by high Sersic 

indexes (n > 4: iKormendv et al' 20091) . Observations of blue 


colour gradien ts from the centres of early-type galaxies towards 
their outskirts (iPeletier et al.lll990t iLiu et al.ll20()5t 


Rudick et al 


201T), mainly attributed to a gradient in me t allictv (iTamura et af. 
200C; iLoubser & Sanchez-Blazaue3l2012^ iMontes et al.l 20l3) 


are also in agreement with a change in stellar properties driven 
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SUBARU Telescope under program SlOA-039. 


by the accretion of smaller systems. Records of accretion events 
are also revealed in the form of spatially extended low-surface 
brigh t ness features in the outsk i rt of stellar halos (iMihos et al.l 
120051: Ivan Dokkum et all 120141: iDuc et al] l2015l) . These sub¬ 
structures are not in a phase-mixed equilibrium in the host 
galaxy potential and, therefore, can in principle be traced as 
kinematic features in velocity phase-space. 


M87 is one of the nearest central galax i es (at an adopted 
distan ce of 14.5 Mpc; ICiardullo et al.l 1200^ iLongobardi et al.l 
|2015[). close to the dy namical centre of the Virgo clus- 
ter ( Binggeli et al .l 119871). It is consid ered a type-cD galaxy 
(I Weil et al.ll997t Kormendv et al.l200^ with an extended stellar 
envelope that reaches ~ 150 kpc in radius. A blue colour gradi¬ 
ent towards its outskirts (iRudick e t alJ201(ih has been interpreted 
as ag e and metallicity gradients dLiu et akl 120051 : iMontes et al.l 
120141) . consistent with a late build-u p of its halo. The or bital 
properties of globular clusters (GCs) (Agnello et aDl2014l) and 
ultra compact dwarfs (IZhang et al.ll2015 1 also favour accretion 
onto the hal o. M87 has been the target of severa l deep imag¬ 
ing surveys (iMihos et al.l2005l:lRudick et alj|2010l) . and its close 
proximity has made it possible to identify hundreds of plane¬ 
tary neb ulas (PNs), and to measure their line-of-sight velocities 
(LOSVs: lDohertv et al.ll200^ ILongobardi et al.ll2015l) . 


In this letter we use the synergy between PN kinematics and 
deep imaging to identify an on-going accretion event in the outer 
halo of M87, which we find to be a non-negligible perturbation 
of the galaxy properties at the distances where it is traced. 
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Fig. 1. Top-left; Projected phase-space, Vlos vs. major axis distance (R), for all spectroscopically confirmed PNs (black asterisks) 
in the halo of M87. Red lines separate the elliptical annuli analysed to isolate cold components associated to the substructure. 
Top-right: Histograms of the LOSVD in the three elliptical annuli. In each panel, the blue lines show the best-fit model computed 
as a combination of three Gaussians. Black-dashed lines show the relative contribution of each component to the LOSVD, with 
parameters for the cold components given in the plot. Bottom-left; As in the top-left panel, however, the green circles and magenta 
diamonds show th e PNs associated with the c old secondary peaks in the LOSVD. Orange squares show a kinematically selected GC 
substructure from iRomanowskv et al.l (l2012h . Bottom-right; Probability that a PN is drawn from the halo component (dark grey 
area) or from the chevron (green, magenta areas). Stars represent PN probabilities at their measured Vlos- 


2. Kinematic evidence for an accretion event in M87 


We acquired kinematic data from the FLAMES/VLT spectro¬ 
scopic survey_for_ajM;ges^pl£_of PNs in the outer regions of 
M87 (iLongobardi et al.ll2013l l2015l) . The total sample consists 
of 254 objects classified as M87 halo PNs and 44 intracluster 
PNs, for which we have obtained LOS Vs with an estimated me¬ 
dian velocity accuracy of 4.2 kms ' . We concentrate here on the 
M87 halo PNs that cover a range of radii from ~ 15 - 150 kpc. 

We see a notable chevron (or “V” shape) structure in the pro¬ 
jected phase-space of the PN sample as shown in Fig. [1] (top 
left). To isolate this kinematical substructure we utilise a three- 
component Gaussian mixture model to identify high-density, 
narrow features on top of a broader distribution. We note that 
there is not enough data to statistically favour this model over 
simpler models with B IG or AIC (Bayesian/Akaike Information 
Criteria: [Liddlell2()(y^ : however, it is visually indicated and we 
will confirm it with photometry in SecH] A brief description of 
the technique is given in t he following paragraph ; for more de¬ 
tails we refer the reader to lPedregosa et al.l (1201 ih . 

A Gaussian Mixture Model (GMM) is a probabilistic model, 
which assumes that a distribution of points can be described as a 


linear combination of K independent Gaussian probability den¬ 
sity functions (PDFs), or components, expressed by: 

K 

p{x) = ^Pk{x\pk,o-k)Pk. (1) 

*=1 

where, v is a data vector (here the LOS Vs), Pk is the mixture 
weight that satisfies the conditions 0 < Pa < 1 and 2f=i 
and p{x\iJ.k,crk) are the individual Gaussian PDFs, with mean 
fill, and dispersion cTk- The GMM classifier implements the 
Expectation-Maximization (EM) algorithm, i.e. an iterative pro¬ 
cess that continuously updates the PDF parameters until conver¬ 
gence is reached. At the end of the EM procedure, the posterior 
probabilities, yk{x) for a data value to belong to each of the k 
Gaussian components are returned. These are described by: 

Pk{x\pk,o-k)Pk ,,, 

Jkix) =-—-. ( 2 ) 

p{x) 

To apply the GMM to our LOSV distribution (LOVSD) we bin 
the PN M87 halo sample into seven elliptical annuli, or stripes 
in phase-space, covering the entire PN velocity phase-space. The 
LOSVD in each annulus is analysed as a combination of three 
Gaussians, where the centres ipk), widths (cta), and weights (Pa) 
are treated as free parameters in the EM algorithm, and have 
uncertainties = ct aV xfS, cTrr,. - cT k! and crp^ — Pk I VS, 
with S - [2„ jkixn)] (lMacKavlf2003h . 
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We find cold components in three out of seven elliptical bins, 
for which we show the histogram of the data, along with the best- 
fit GMM and reduced in Fig. [T| (top-right panel). We also 
plot (bottom-left panel) the LOSV phase-space for the 254 PNs 
in the halo of M87: black crosses represent PNs of the smooth 
halo LOSVD of M87, while magenta diamonds and green dots 
are PNs that have a higher probability (see eq. |2l) to belong to 
the chevron. Finally, we show the probability that a given PN 
is drawn from each of the components as a function of its ve¬ 
locity (bottom-right panel). The GMM assigned a total of 54 
PNs to the chevron substructure, which covers 700" (~ 50 kpc) 
for major axis distances 500" < R < 1200". The separation, 
AV, between the two peaks of the cold components becomes 
smaller at larger distances. For the three elliptical bins, it is 
AV = 654.5 + 40.1, 246.4 + 20.5, and 153.3 + 20.6kms-', 
respectively. At R ~ 1200" (~ 90 kpc) the width of the chevron 
goes to zero with L OS Vs close to the galaxy ’s systemic velocity 
(Vsys = 1275 kms~*: iLongobardi et all 2015[n . PNs on the arms 
of the chevron are seen on both the northern and southern sides 
of the galaxy as is shown in Fig. (see Sect. O. The broad 
Gaussian with average mean velocity ~ 1290 kms“^ and disper¬ 
sion ~ 320 kms ' in the three bins traces the M87 halo. 

The search for kinematic features in the phase-space of GCs 
has resulted in the discovery of a similar chevron-like struc¬ 
ture (iRomanowskv et aT]l2012h . shown in Fig[T] above as orange 
squares. Though the morphology in the phase-space is similar 
it differs in a number of physical properties; the width goes 
to zero at Rqc ~ 1500" with Vlosgc = 1307 kms^', versus 
R ~ 1200" and Vlos ~ 1250 +21 kms-' for the PNs. Moreover, 
the 27 chevron GCs show a very different spatial distribution 
with the highest density of points on the NE photome tric minor 
axis (IRomanowskv et al1l2012t ID’ Abrusco et al.ll2015h . and few 
GCs near the crown substructure traced by the PNs. 

3. Localising the substructure with deep imaging 

In Fig. 12 (left), we show the position of the chevron PNs over¬ 
plotted on 1.6x 1.6 deg“ V-band image of M87, with an estimated 
surfac e brightness limit of /ry = 28.5 magarcsec-^ (iMihos et alJ 
l2015l) . Here we see the large spatial extent of the substructure 
associated with the chevron PNs. Because of its shape on the 
image we will refer to it as the veil of M87. Now we are inter¬ 
ested to see if this feature is also visible in the optical light. 

To this end, we constructed an unsharped masked im¬ 
age, which is the difference between the original image and a 
smoothed image. We utilised the IRAF task fmedian to smooth 
the original image by using a window with a size of 1450" x 
1450". This window size was chosen so that it contained the 
large scale extension of the substructure. By looking at the high¬ 
est concentration of PNs, in the NW region of M87, it can be 
seen to extend over many hundreds of arcseconds (~ 800"). 
Thus, the adopted box size is ~ 1.8 times larger than the long 
side of the feature. 

The results of the unsharped masking can be seen in Fig. |2 
(central panel), where the high frequency structures are now 
clearly visible. A previously unknown debris structure, with a 
crown-like shape, can be seen on top of M87 at the NW side. 
This we refer to as the crown of M87’s veil. It has a character¬ 
istic width of ~ 300", an extension of ~ 800", and is almost 
perpendicular to M87’s photometric major axis. Just as the PN 

' The stability of the fitted parameters and the measured distance of 
the chevron edge were tested with 100 GMM runs for different mock 
data sets and initialisation values. 


spatial distribution showed, the edge of this feature is found at 
R~ 1200" (~ 90 kpc). 

Over the same major axis distances at which the substructure 
is located, we also observe variations in the M87 ellipticity pro¬ 
file. Between 300" < R < 800" the ellipticity i ncreases and then 
flattens to a value of e ~ 0.43 for R > 800" (iKormendv et al.l 
l2009l) . In particular, the region at which the gradient flattens re¬ 
flects the crown-like overdensity. 

4. Physical properties of the accreted satellite 

4.1. Luminosity and a-parameter 

To understand the physical origin of this structure, we compute 
its total luminosity. We do it in a region with size ~ 800" x 300", 
after the subtraction of the local background, dete rmined using 
a regi on photometry method (for more details see iRudick et al.l 
I 2 OIOI) . We find a total luminosity in the range of Laown = 3.7 + 
0.9 X IO^Lq.v- When compared to the luminosity determined in 
the same region from the original image, L ~ 5.9 x IO^Lq v, 
it is clear that at these distances the substructure represents a 
significant fraction, ~ 60%, of the total light. 

In the region of the crown we count Npn, crown = 12 + 3 PNs, 
while we find a total of NpN.chevron = 54 + 7 PNs associated to 
the entire chevron (see Sect. O . By correcting these nu mbers 
for incompleteness factors as in iLongobardi et al.l (12015h these 
become 19 and 142, respectively. Hence, by scaling Lcrown to 
the total number of PNs associated to the chevron we obtain the 
total luminosity associated to the progenitor of the M87 veil to 
beL~2.8 + 1.0xl0®Lo,v. 

The total number of PNs is proportional to the total bolo- 
metric luminosity of the parent stellar population, and the pro¬ 
portionality is quantifi ed with the luminos ity-specific PN den¬ 
sity, or, cr-parameter (iBuzzoni et al.ll200^ . Utilising the com¬ 
puted luminosity, and the completeness-corrected estimate for 
NpN,chevron, we Can Calculate the cr-parameter for the progenitor 
of the substructure. Considering that the typical probability of 
the NpN,chevron PN to belong to the chevron is ~ 0.7 (Fig.[TJ, we 
obtain Q' 2.5 = 1.8 + 0.7 x IQ-^NpnL -'^^,. Here 02.5 is 2.5 mag 
down the luminosity function as in ILongobardi et ^ (12015l) . 
and we have assumed a bolometr ic correction for the V-band 
of BCv=0.85 (|Buzzoni_et_^ ]l2006ll and BCo=-0.07 for the Sun. 

4.2. Coiour and Mass 

In Fig. 12 (right panel) we show the chevron PNs overplotted on 
the B-V colour image of M87, that combines the V-band data 
(see Sect.O with deep B-band photometry with a surface bright¬ 
ness limit of pb = 29 mag arcsec-^. It is interesting to notice that 
close to R ~ 1200" the colour shows an azimuthal variance, 
such that alon g the photometric m inor axis the measured val¬ 
ues are redder (iMihos et alJl20L5h . This feature correlates with 
the spatial number density of the chevron PNs, showing a deficit 
in number along the photometric minor axis. This suggests that 
the bluer regions are the result of the accreted material on top 
the light from M87’s halo. In particular, the crown structure is 
measured to have integrated colour (B - V) = 0.76 + 0.05. 

(B-V) colour is a good estimator of the mass-to-light- 
ratio, T*, of the underlyi ng stellar population. By ado pting 
T* = 2.3 for (B-V)=0.76 (iMcGaugh & Schornbertllml . the 
total stellar mass associated to the disrupted galaxy is then 
M = 6.4 + 2.3 X lO^Mo. 

From the distribution and velocities of chevron PNs in Fig.|2 
a possible interpretation of the satellite orbit could be that it was 
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Fig. 2. Spatial and colour distribution associated with the kinematic su bstructure id e ntified in the phase-space of the M87 halo PNs. 
Left panel: V-band image of a 1.6 x 1.6 deg^ centred on M87 from iMihos et al.l (l2015h . Full circles, and diamonds indicate the 
spatial position of the M87 halo PNs in the chevron substructure. Magenta and green colours indicate PN LOSVs above and below 
Vlos = 1254kms the LOSV at the end of the chevron. Central panel: Unsharp masked image of M87 median binned to enhance 
faint structures. The crown-shaped substructure is visible at distance of 800"-1200"(~ 60 - 90 kpc) along the major axis, NW of 
M87. Details are given in Section[3] The blue line measures 90 kpc. Right panel: (B-V) colour image of M87 from iMihos et al.l 
(1201 5h with chevron PNs overplotted (white dots). The dashed ellipse indicates the isophote at a major axis distance of 1200". The 
crown is found in a region where the (B-V) colour is on average 0.8, bluer than on the minor axis. 


first disrupted entering M87 from the South (along the green 
dots), with the debris then moving up North, turning around in 
the crown region, and coming back down on both sides across 
M87 (the veil, magenta dots). The velocities would then imply 
that the northern side of M87 is closer to the observ er. The dy¬ 
namic al time for such an orbit is of order < 1 Gyr (IWeil et alJ 
l1997h . 


5. Summary and Conclusion 

In this letter we have presented kinematic and photometric evi¬ 
dence for an accretion event in the halo of the cD galaxy M87. 
This event is traced by PNs whose velocity phase-space shows 
a distinct chevron-like feature, which is a result of the incom¬ 
plete phase-space mixing of a disrupted galaxy. At major axis 
distances ofR~ 60-90 kpc, where the width of the chevron goes 
to zero, a deep optical image shows the presence of a crown-like 
substructure that contributes > 60% of the total light in this area. 

The crown of M87’s veil is the densest part of the entire sub¬ 
structure, which covers ~ 50 kpc along the major axis. In this 
region also a radial variation in M87’s ellipticity profile is ob¬ 
served. Looking at the spatial distribution of all the chevron PNs, 
it traces the azimuthal variation observed in the colour of M87, 
showing a deficit in number of tracers along the photometric mi¬ 
nor axis where the galaxy is redder, and a higher fraction where 
the substructure is strongest and the colour is bluer. 

We determined several physical properties of the disrupted 
satellite: a total luminosity of L = 2.8 + 1.0 x IO^Lq v, 
colour (B-V)=0.76 + 0.05, and total stellar mass of 
M = 6.4 + 2.3 X 10^ Mq. The inferred value for the a-parameter 
is O' = 1.8 + 0.7 X 10“® NpNLg'j^^j. The similar colours of the 
accreted satellite and ICL suggest that the cD halo of M87 
is presently growing by the accretion of similar star-forming 
systems as those that originate the diffuse IC component. 

The evidence for on-going accretion in the outer halo of M87 
is consistent with the observed size growth of giant elliptical 
galaxies and with predictions by theory. The presence of the 


newly discovered substructure within the halo of M87 demon¬ 
strates, that beyond a distance of ~ 60 kpc, its halo is still as¬ 
sembling. 

Acknowledgements. AL is grateful to J. Elliott for helpful comments which im¬ 
proved the manuscript. JCM is supported by NSF grant 1108964. 
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